#File Name: FullWDITransparencyIndexConstruct.r
#Author: James Hollyer
#Date: 03/03/2013
#Purpose: To construct a measure of government transparency.  To do this, we
# rely on the reporting or missingness of data from the WDI index.  
#Data Input: IRTprep.dta
#OS: Windows 7

#The following will load the data into R.

library(R2jags)
library(rjags)
library(foreign)

setwd("c:/users/james/desktop/dropbox/transparency_and_democracy/transparencyindex/PAReplicationMaterials/PostPred Checks/")
WDIdata<-read.dta("c:/users/james/desktop/dropbox/transparency_and_democracy/transparencyindex/PAReplicationMaterials/FullWDIIRTprep2.dta")

attach(WDIdata)

#The following will prepare a data matrix for export to BUGS.


data.b<-list("var_1", "var_2", "var_3", "var_4", "var_5",
"var_6", "var_7", "var_8", "var_9", "var_10", "var_11", "var_12", "var_13",
"var_14", "var_15", "var_16", "var_17", "var_18", "var_19", "var_20", "var_21",
"var_22", "var_23", "var_24", "var_25", "var_26", "var_27", "var_28", "var_29",
"var_30", "var_31", "var_32", "var_33", "var_34", "var_35", "var_36", "var_37",
"var_38", "var_39", "var_40", "var_41", "var_42", "var_43", "var_44", "var_45",
"var_46", "var_47", "var_48", "var_49", "var_50", "var_51", "var_52", "var_53",
"var_54", "var_55", "var_56", "var_57", "var_58", "var_59", "var_60", "var_61",
"var_62", "var_63", "var_64", "var_65", "var_66", "var_67", "var_68", "var_69",
"var_70", "var_71", "var_72", "var_73", "var_74", "var_75", "var_76", "var_77",
"var_78", "var_79", "var_80", "var_81", "var_82", "var_83", "var_84", "var_85",
"var_86", "var_87", "var_88", "var_89", "var_90", "var_91", "var_92", "var_93",
"var_94", "var_95", "var_96", "var_97", "var_98", "var_99", "var_100", "var_101",
"var_102", "var_103", "var_104", "var_105", "var_106", "var_107", "var_108",
"var_109", "var_110", "var_111", "var_112", "var_113", "var_114", "var_115",
"var_116", "var_117", "var_118", "var_119", "var_120", "var_121", "var_122",
"var_123", "var_124", "var_125", "var_126", "var_127", "var_128", "var_129",
"var_130", "var_131", "var_132", "var_133", "var_134", "var_135", "var_136",
"var_137", "var_138", "var_139", "var_140", "var_141", "var_142", "var_143",
"var_144", "var_145", "var_146", "var_147", "var_148", "var_149", "var_150",
"var_151", "var_152", "var_153", "var_154", "var_155", "var_156", "var_157",
"var_158", "var_159", "var_160", "var_161", "var_162", "var_163", "var_164",
"var_165", "var_166", "var_167", "var_168", "var_169", "var_170", "var_171",
"var_172", "var_173", "var_174", "var_175", "var_176", "var_177", "var_178",
"var_179", "var_180", "var_181", "var_182", "var_183", "var_184", "var_185",
"var_186", "var_187", "var_188", "var_189", "var_190", "var_191", "var_192",
"var_193", "var_194", "var_195", "var_196", "var_197", "var_198", "var_199",
"var_200", "var_201", "var_202", "var_203", "var_204", "var_205", "var_206",
"var_207", "var_208", "var_209", "var_210", "var_211", "var_212", "var_213",
"var_214", "var_215", "var_216", "var_217", "var_218", "var_219", "var_220",
"var_221", "var_222", "var_223", "var_224", "var_225", "var_226", "var_227",
"var_228", "var_229", "var_230", "var_231", "var_232", "var_233", "var_234",
"var_235", "var_236", "var_237", "var_238", "var_239", "var_240")


#List the parameters to be estimated

parameters<-c("beta1", "beta2", "beta3", "beta4", "beta5", "beta6", "beta7", "beta8", "beta9",
"beta10", "beta11", "beta12", "beta13", "beta14", "beta15", "beta16", "beta17",
"beta18", "beta19", "beta20", "beta21", "beta22", "beta23", "beta24", "beta25",
"beta26", "beta27", "beta28", "beta29", "beta30", "beta31", "beta32", "beta33",
"beta34", "beta35", "beta36", "beta37", "beta38", "beta39", "beta40", "beta41",
"beta42", "beta43", "beta44", "beta45", "beta46", "beta47", "beta48", "beta49",
"beta50", "beta51", "beta52", "beta53", "beta54", "beta55", "beta56", "beta57",
"beta58", "beta59", "beta60", "beta61", "beta62", "beta63", "beta64", "beta65",
"beta66", "beta67", "beta68", "beta69", "beta70", "beta71", "beta72", "beta73",
"beta74", "beta75", "beta76", "beta77", "beta78", "beta79", "beta80", "beta81",
"beta82", "beta83", "beta84", "beta85", "beta86", "beta87", "beta88", "beta89",
"beta90", "beta91", "beta92", "beta93", "beta94", "beta95", "beta96", "beta97",
"beta98", "beta99", "beta100", "beta101", "beta102", "beta103", "beta104",
"beta105", "beta106", "beta107", "beta108", "beta109", "beta110", "beta111",
"beta112", "beta113", "beta114", "beta115", "beta116", "beta117", "beta118",
"beta119", "beta120", "beta121", "beta122", "beta123", "beta124", "beta125",
"beta126", "beta127", "beta128", "beta129", "beta130", "beta131", "beta132",
"beta133", "beta134", "beta135", "beta136", "beta137", "beta138", "beta139",
"beta140", "beta141", "beta142", "beta143", "beta144", "beta145", "beta146",
"beta147", "beta148", "beta149", "beta150", "beta151", "beta152", "beta153",
"beta154", "beta155", "beta156", "beta157", "beta158", "beta159", "beta160",
"beta161", "beta162", "beta163", "beta164", "beta165", "beta166", "beta167",
"beta168", "beta169", "beta170", "beta171", "beta172", "beta173", "beta174",
"beta175", "beta176", "beta177", "beta178", "beta179", "beta180", "beta181",
"beta182", "beta183", "beta184", "beta185", "beta186", "beta187", "beta188",
"beta189", "beta190", "beta191", "beta192", "beta193", "beta194", "beta195",
"beta196", "beta197", "beta198", "beta199", "beta200", "beta201", "beta202",
"beta203", "beta204", "beta205", "beta206", "beta207", "beta208", "beta209",
"beta210", "beta211", "beta212", "beta213", "beta214", "beta215", "beta216",
"beta217", "beta218", "beta219", "beta220", "beta221", "beta222", "beta223",
"beta224", "beta225", "beta226", "beta227", "beta228", "beta229", "beta230",
"beta231", "beta232", "beta233", "beta234", "beta235", "beta236", "beta237",
"beta238", "beta239", "beta240")

#send to bugs
results <- jags(data=data.b, inits=NULL, parameters, model.file="ConstantsOnly.txt",
n.chains=2, n.iter=10000, n.burnin=5000)

print(results)

attach.jags(results)
TePCP<-cbind(1:500)

for (j in 1:240){
   eval(parse(text=sprintf("diff%s <- cbind(1:3875)", j)))
}

for (t in 1:500){
TePCP[t] <- 0
  for (j in 1:240){
  eval(parse(text=sprintf("p.hat%s <- 1/(1+exp(-beta%s[%s]))",j,j,t)))
    for (c in 1:3875){
    eval(parse(text=sprintf("diff%s[%s] <- 1 - abs(WDIdata$var_%s[%s] - p.hat%s)",j,c,j,c,j)))
    }
  eval(parse(text=sprintf("TePCP[t] <- TePCP[t] + mean(diff%s)", j)))
  }
}
save.image("ConstantsOnly.RData")





